In [ ]:
import kwant
import tbmodels
import wraparound
import numpy as np
import scipy.linalg as la
import scipy.sparse.linalg as sla
import matplotlib
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
from matplotlib import pyplot
import csv
from mayavi.mlab import *
from mayavi import mlab
import mayavi
import math
Sample model for number of neighbor checks (nn)¶
In [2]:
model = tbmodels.Model.from_wannier_files(
hr_file='BiPentaSi_hr.dat',
wsvec_file='BiPentaSi_wsvec.dat',
xyz_file='BiPentaSi_centres.xyz',
win_file='BiPentaSi.win'
)
nn = 2
remove = [R for R in model.hop if (R[0] > nn) or (R[1] > nn) or (R[2] > nn)]
for R in remove:
del model.hop[R]
remove = [R for R in model.hop if (R[0] < -nn) or (R[1] < -nn) or (R[2] < -nn)]
for R in remove:
del model.hop[R]
lattice = model.to_kwant_lattice()
def lead_shape(p):
x, y, z = p
return -1*model.uc[0][0] < x < 1*model.uc[0][0] and -1*model.uc[1][1] < y < 1*model.uc[1][1] and 0 < z < 30
sym_lead = kwant.TranslationalSymmetry(
lattice.vec((1, 0, 0)),
lattice.vec((0, 1, 0)),
lattice.vec((0, 0, 1))
)
lead = kwant.Builder(sym_lead)
lead[lattice.shape(lead_shape, (0, 0, 15))] = 0
model.add_hoppings_kwant(lead)
kwant_sys = wraparound.wraparound(lead).finalized()
fig = kwant.plot(kwant_sys, site_size=0.18, site_lw=0.01, hop_lw=0.05, site_color='g', dpi = 300)
Dispersion curve for 1x1 unit cell in the kx direction¶
In [3]:
k_list = [(kx, 0, 0) for kx in np.linspace(0, 0.5, 50)]
kx = range(50)
eigs_kwant = [la.eigvalsh(
kwant_sys.hamiltonian_submatrix(2 * np.pi * np.array(k))
) for k in k_list]
fig, ax = plt.subplots()
for band in np.array(eigs_kwant).T:
ax.plot(kx, band+1.784, 'b')
plt.ylim(-3,3)
Out[3]:
6x2 model for transport calculations¶
In [4]:
def shape(p):
x, y, z = p
return -3*model.uc[0][0] < x < 3*model.uc[0][0] and -1*model.uc[1][1] < y < 1*model.uc[1][1] and 0 < z < 30
def lead_shape(p):
x, y, z = p
return -3*model.uc[0][0] < x < 1*model.uc[0][0] and -1*model.uc[1][1] < y < 1*model.uc[1][1] and 0 < z < 30
remove = [R for R in model.hop if (R[0] > nn) or (R[1] > nn) or (R[2] > nn)]
for R in remove:
del model.hop[R]
remove = [R for R in model.hop if (R[0] < -nn) or (R[1] < -nn) or (R[2] < -nn)]
for R in remove:
del model.hop[R]
lattice = model.to_kwant_lattice()
kwant_model = kwant.Builder()
kwant_model[lattice.shape(shape, (0, 0, 15))] = 0
model.add_hoppings_kwant(kwant_model)
print("Model for nn=",nn,"created")
#1D translational symmetry
sym_lead = kwant.TranslationalSymmetry(
lattice.vec((-3, 0, 0))
)
lead = kwant.Builder(sym_lead)
lead[lattice.shape(lead_shape, (0, 0, 15))] = 0
model.add_hoppings_kwant(lead)
print("Lead for nn=",nn,"created")
#Add the norbs in the system for plotting of local quantities (current, electron density)
for x in kwant_model.sites():
x.family.norbs=1
for x in lead.sites():
x.family.norbs=1
kwant_model.attach_lead(lead)
kwant_model.attach_lead(lead.reversed())
kwant_sys = kwant_model.finalized()
kwant.plot(kwant_sys, site_size=0.5, site_color='g', file='wire-Xdirection.png', dpi = 600)
Out[4]:
In [5]:
Ef = -1.784
energies = [i+(Ef-1) for i in np.linspace(0,2,100)]
data = []
for energy in energies:
print(energy)
# compute the scattering matrix at a given energy
try:
smatrix = kwant.smatrix(kwant_sys, energy)
data.append(smatrix.transmission(1, 0))
except RuntimeError:
print("Error")
continue
In [6]:
pyplot.figure()
energies = np.array(energies)
pyplot.plot(energies-Ef, data)
pyplot.xlabel("E-Ef (eV)")
pyplot.ylabel("Conductance (2e^2/h)")
pyplot.show()
2D current plot at the z location of Si7 and Si8 (zcut=35) at the energy where their ldos is at its maximum¶
In [7]:
for x in kwant_model.sites():
x.family.norbs=1
for x in lead.sites():
x.family.norbs=1
kwant_model.attach_lead(lead)
kwant_model.attach_lead(lead.reversed())
kwant_sys = kwant_model.finalized()
energy = Ef-0.191919
#-1.97591919
psi = kwant.wave_function(kwant_sys, energy=energy)
J = kwant.operator.Current(kwant_sys)
all_states = np.vstack(psi(0))
current = sum(J(p) for p in all_states)
curr_i = kwant.plotter.interpolate_current(kwant_sys,current)
In [9]:
for zcut in range(1,45):
print(zcut)
kwant.plotter.streamplot(curr_i[0][:,:,zcut,:2], (curr_i[1][0],curr_i[1][1]),dpi=100)
